
set.seed(7654)

num_obs <- 120
alpha <- 0.05
sims <- 1000
sim <- NULL

possible.taus <- seq(-0.05, 0.05, by = 0.002)


for (j in 1:length(possible.taus)){ 
  
  tau <- possible.taus[j] # Pick the jth value for N 
  
  Y0 <- rnorm(n=num_obs, mean=0.43, sd=0.023) #  potential outcome 
  Y1 <- Y0 + tau # underperformance potential outcome                                   
  significant.experiments <- rep(NA, sims) # Empty object to count significant experiments
  
  #### Inner loop to conduct experiments "sims" times over for each N #### 
  for (i in 1:sims){
    Z.sim <- rbinom(n=num_obs, size=1, prob= 0.5) # Do a random assignment 
    Y.sim <- Y1*Z.sim + Y0*(1-Z.sim) # Reveal outcomes according to assignment 
    fit.sim <- lm(Y.sim ~ Z.sim) # Do analysis (Simple regression) 
    p.value <- summary(fit.sim)$coefficients[2,4] # Extract p-values 
    significant.experiments[i] <- (p.value <= alpha) # Determine significance according to p <= 0.05
  }
  if(!exists("final_df")){
    final_df <- cbind(tau, mean(significant.experiments)) %>% data.frame()
  } else{
    final_df <- bind_rows(final_df, cbind(tau, mean(significant.experiments)) %>% data.frame())
  }
}

main_power <-
  final_df %>%
  ggplot(aes(x = tau*100, y = V2*100)) +
  geom_point(size = 1) +
  geom_line() +
  #geom_vline(xintercept = 3000, color = "red", linetype = "dashed") +
  geom_vline(xintercept = 0, color = "pink", linetype = "dashed") +
  geom_hline(yintercept = 80, color = "pink", linetype = "dashed") +
  #geom_hline(yintercept = 80, color = "red", linetype = "dashed") +
  theme_bw() +
  ylab("Probability of Detecting Effect (%)") +
  xlab("Treatment Effect\n(Share Resources, Percentage Points)") +
  theme(panel.grid.minor = element_blank(), 
        panel.grid.major.x = element_blank(), 
        panel.grid.major.y = element_line(size = .2)) + 
  theme(strip.background =element_rect(fill="white")) +
  theme(panel.border = element_blank(),
        axis.line.x = element_line(size = 0.5, linetype = "solid", colour = "black"),
        axis.line.y = element_line(size = 0.5, linetype = "solid", colour = "black"))


ggsave(plot = main_power, "./outputs/figures/figure_a3.pdf", width = 5, height = 3)


